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A NUMERICAL METHOD FOR DETERMINING THE NATURAL 


VIHRATION CHARACTERISTICS OF ROTATING NONUNIFORM CANTU EVER BLADES 

William F. White, Jr. 

Langley Directorate 
U.S. Army Air Mobility R&D Laboratory 
NASA Langley Researi'h Center 
Hampton, Virginia 23665 

and 

Raymond E. Malatino 
Naval Air Systems Command 
V/ashington, D.C. 20361 

SUMMARY 

A method is presented for determining the free vibration characteristics 
of a rotating blade having nonuniform spanwise properties and cantilever 
boundary conditions. The equations which govern the coupled llapwise, chord- 
wise, and toi^sional motion of such a blade are solved using an integrating 
matrix met'..._. By expressing the equations of motion in matrix notation, 
utilizing the integrating matrix as an operator, and applying the boundary 
conditions, the equations are formulated into an eigenvalue problem whose 
solutions may be determined by conventional methods. Computed results are 
compared with experimental data. 

INTRODUCTION 

' Natural vibration characteristics of rotating blades are of fundamental 
importance from the viewpoint of flying qualities, blade life, vibration 
levels, and stability. Helicopters, propellers, and turbines may have serious 
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resonant vibration problems when the excitation frequencies are equal to some 
multiple of the rotational speed. To insure that conditions susceptible to 
resonance do not exist within the range of operating speeds, it is necessary 
that the natural frequencies be determined accurately. Also, the natural 
modes, because of their orthogonality relationships, are often used in forced 
response and stability calculations. 

This paper formulates a numerical solution of the natural vibration fre- 
quencies and mode shapes of rotating nonuniform blades. This problem has 
been treated analytically in a very complete development by Houbolt eind 

Brooks.^ However, very few results are presented, and they are for special 

2-lU 

cases of limited interest. Numerous other studies have investigated 
various facets of the problem using a variety of nvunericed solution methods. 
Generally, these studies may be classified as investigations of a subset of 
the governing equations derived in Reference 1. The present analysis employs 
the governing equations derived in Reference 1. The integrating matrix as 
developed by Hunter^^ is the basis for the present method of solution. The 
integrating matrix is a means of numerically integrating a function that is 
expressed in terms of the values of the function at inc ements of the inde- 
pendent variable. It is derived by expressing the j.ntegrand as a polynomial 
in the form of Newton’s forward-difference interpolation formula. 

The equations of motion which are linear homogeneous equations having 
variable coefficients ere expressed in matrix form using the integrating 
matrix. The boundary conditions are applied and the resulting matrix equa- 
tions expressed in standard eigenvalue form. Solutions to this eigenvalue 
problem may be obtained by conventional methods. In developing the solution, 
all functions are in effect represented by seventh-degree polynomials at the 
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boundaries as well as elsewhere on the beam. Since the polynomials approxi- 
mate the functions very accurately, the integration of these polynomial 
representations yield extremely small errors. The method is appealing be- 
cause the numerical solution may be formulated quickly from the governing 
equations and may be easily programmed for computations by a digital computer. 
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SYMBOLS 


Blade section constants 

Young's modulus of elasticity 
Ratio of tip to root bending modulus 

Distance between mass and elastic axis , positive when mass 
axis lies ahead of elastic axis 

Distance between area centroid of tensile member and elastic 
axis, positive for centroid forward 

Shear modulus of elasticity 

Cross-section area moments of inertia 


Torsional stiffness constsmt 


Polar radius of gyration of cross-sectional area effective . 1 
carrying tension 


Polar radius of gyration of cross-sectional mass about elastic 


axis 


(k. 


2 = k 


m 


m. 


k 


m. 


Mass radii of gyration about major neutral axis and about an 
axis perpendicular to chord through the elastic axis, 
respectively 

Resultant moments in y - and z - directions, respectively 


Mass per unit length 
Ratio of tip to root mass 
Number of blade stations 

Resultant cross-sectional torque about elastic axio 
Blade radius 
Blade tension, T = 



mx dx 
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i V ,v 

y z 

v,w 

Cross-sectional shears in y- and z-directions 

Lateral displacements of beam, in plane of rotation and normal 
to plane, respectively 

v,w 

Vibration amplitude of v and w, respectively 

x,y ,2 

1 

« 

Coordinate system which rotates with blade such that x-axis 
falls along initial or undeformed position of elastic axis 

n 

Radial position at which applied loadings are assumed to act 

e 

Blade angle prior ■' any deformation, positive when leading 
edge is up 

* X 

Eigenvalue 

$ 

Total torsional deflection, 4* = eigenvector 


Elastic torsional deflection, positive when leading edge is up 

<t> 

Steady-state twist (Appendix B) 

i 4) 

Vibration amplitude of <j> 

! n 

1 

Rotor angular velocity 

0) 

Natural frequency of vibration 

Primes denote derivatives with respect to x; dots denote derivatives with 
respect to time 

Matrix 

notation: 

[ ] 

Square matrix 

r j 

Diagonal matrix 

{ } 

Column matrix 

: [ 

Inverted matrix 

• . [1] 

Unit or identity matrix 
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EQUATIONS OF MOTION 

The beam considered in the present analysis is shown in Figure 1. This 
case represents the coupled bending and torsion of a twisted rotating beam where 
the elastic axis, mass axis, and tension axis are not necessarily coincident. 

The free vibration equations of motion for combined bending and torsion are 
derived in Appendix A in tVe form 



cos 


e + Elg 



sin^ 0)w" (Elg 
♦ cos 6 + 

(mw + me$ cos 6) 



- EI^) cos 6 sin 9 v" - EB20'4»* sin 0 
{mn[w(n) - w(x)] + men cos 0} dn 
(n - x)dn 


(la) 


(El, - 


EI^) cos 0 sin 0 w" + (EI^ sin^ 0 + Elg cos^ 9)v" - EBg 0'i^' cos 9 




+ Te (j) sin 0 + fl I {m[xv(n) - nv(x)] - mencfi sin 0} dn 

& _/ X 




(mv - me$ sin 0) (n - x)dn 


(lb) 


[GJ + Tk^^ + EB^ (0')^]‘t>' - EBg 0' (v" cos 0 + w" sin Q) + f Xe (v" sin 0 

•'x ^ 

2 /*^ 2 2 
- w" cos 0) dn + / [mev sin 0 + m (k - k )d> cos P0 

Jx “2 “l 

= - [me(w 


- men (V sin 0 - w' cos 0)] dn 


cos 0 - V sin 0 ) 


+ 4>m k ] dn 
n 


(Ic) 
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where the variation of the axial tensile force Is given by 

2 

T' + m n X ■ 0 ^2 

Also, V, V, and 4> o.re the edgewise, flapwise, and torsional displacements, 
respectively. The principal assumptions used in deriving this system of 
equations are briefly outlined: 

(a) The beam is attached to a rigid hub and rotating at constant speed. 

(b) The beam elasticity is adequately described by the conventional 
bending and torsion oharacterlatics described in Reference 1. Furthermore, 
shear deformation and rotary inertia are assumed to be negligible. 

(c) The elastic axis of the undeformcd bean is a straight line. 

(d) The pitch axis is coincident with the elastic axis of the undefonned 


bean. 

(e) The beam is assumed to have zero precone and prelag. 

(f) The cross section is symmetric about the major principal axis. 

(g) Cyclic pitch is negligible. 

For the present analysis, cantilever- free boundary conditions are assumed. 

Thus, displacements and slopes are zero at the root: 

v(0,t) = w(0,t) * <|)(0,t) = v'(0,t) =» w'(0,t) = 4>'(0,t) = 0 (3) 

also, moments and shears are zero at the tip: 

Q(R,t) = M^(R,t' = M (R,t) = V (R,t) = V (R,t) = 0 (U) 

j z y z 

The integrals in Equation (l) may be evaluated by using a matrix operator [L] 
defined as 



f(x) dx 


= [L]{f} 


(5) 
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Thus, the premixltlplication by [L] of a colunn matrix of the function f(x) 
yields the integration of f(x) from x^ to x^. The (n ♦ l) matrix [l 1 fs 
given by Hunter^^ for polynomials of degrees one to seven. For the present 
analysis, a seventh degree polynomial approximation is used. The matrix 
iLj avoids solving a set of clmultaneous equation? since it Jn<rllcitly contains 
the polynomial coefficients at each station. By expressing the equations of 
motion in matrix notation, utilizing the integrating matrix as an operator, and 
applying the boundary conditions, the equations are formulated into an eigen- 
value problem whose solutions may be determined by conventional methods. 

The following relationships are valid at each spanwise station: 


v*(w,t) « v'(0,t) 

* / v"(x,t) dx; v(x,t) ■ v(0,t) ♦ f v'(0,t) dx 

(6a) 


•'0 Jo 





w'(x,t) ■ w'(0,t) 

J 0 w(x,t) ■ w(0,t) +J w'(w,t) dx 

(6b) 

(j)(x,t) = ((»(0,t) + 

(j)'(x,t) dx 

(6c) 


Applying the boundary conditions given by Equation (3) to Equation (6) and 
using Equation (5) yields 

{v'} - [F]{v"}, {v} = [F] {v*} = [F]2 {v"} 

{w’> = [F]{w"}; {w} = [F]{w*} = [F]^ {w" } 

{(J)} = [F]{<|)'} 

where [F] « til [L] til. The first element in the above column vectors is the 
corresponding quantity at the tip. Assvunlng simple harmonic motion 
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w(x,t) ■ w(x)e^‘**^, v(x,t) ■ v(x)e^'*^, <>{x,t) «« <*’(x)e^‘^ 


and applying the matrix operator [L] to Equation (1) yields 


and 




(7a) 


where 


{♦} 


4 >’ 
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Equation (7a) may be expressed in standard eiger 'alue form as 


X{4>} = [D] {4>} 


(7b) 


Where [D] « [G]~^ [h]. Solutions of Equation (7b) define the natural fre- 


quencies and associated modal vectors. The first element in the above matrices 
corresponds to a tip vedue. In order to satisfy the tip boundary conditions. 
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the 1, n ♦ 1, and 2n + 1 rows and columns of the dynamic matrix [D] are de- 
leted. The partitions of [O] and [h] are (n ♦ 1) matrices defined as follows: 
[Gj ] • tEI rce^ Q + El sin^ 0j + [P (mx)[F]^ 

[G^ ] ■ t(EI^ - El, ) cos 0 sin 0j 
"w 2 ^ 

[G ] " - J^^[L]fe cos 0jtP.(mx)J > n^[L]fmex cos 0j[P] - 9' sin 0j 

3w a 2 2 

[G^ ] « ^{EI 2 - EI^) cos 0 sin 0j 

[Gg ] = fEI^ sin^ 0 + EI^ cos^ 0j - [P^][m] [F]^ + fl^LP^Cmx)] [F]^ 

[G_ ] ■ il^[L]fe sin 0jfP (mx)J + [L]fme sin 0j[F]^ - fi^[L]fmex sin 0J[F] 
a 2 

- fEB^ 0' cos 0j 

C. 

[G- ] * - tP_(mx)Jie cos 0J[F] + ^2^ [L]fmex cos 0j[F] - fEB. 0' sin 0j 

1.2 a 2 

9 

[G_ ] = fP_(mx)Jfe sin 0J[F] - fxJ[L]fme sin 6j[F] - fEB. O' cos Oj 

2.2 a 2 

9 

[G, ] = fGJj + tEB (0')^J + fk ^J?P-(mx)j 4 [L]fm(k^ - k^ ) cos 20j[F] 
3^ 1 a 2 ®2 “l 

[H^ ] = [PjfmJ[F]2 
■‘‘w 

[H ] = [0] 

w 

[H ] = [L]fme cos 0J[F]^ 

[H^ ] = [0] 

V 

[Hg ] = [P^][m][F]^ 

• » 


a 



I 



il 

[n ] * “ [L][ne sin 6][F]^ 

■^v 

[H^ ] = [Pj^] fme cos 6][Fj 

[H^ ] = - [P^)Pme Kin 6j[F] 

4> 

tlL ] - iLltefOlFl 

Where 

[Pj_] » [LltxJ - Px>][L] 

= dii^ [L]{f} 

[P^Cf)] = [L]ffJ - rPgCf)^ 

[F] - U1[L]W 


The m'jaal deflections are determined from the eigenvectors of Equation (7’t)) as 



' Solutions of Equation (7b) are obtained by using the QR transformation 

I 

^ method. Tne dynamic matrix [D] is reduced to upper Hessenberg form by ele- 

^ 15 l6 

mentary similarity transformations. ’ The similarity transformations, 

1 5 IT l6 

known as QP transformations, * * of Francis are used iteratively 

1 to reduce the matrix to an upper triangular form. The eignevectors corresponding 

I 

I to the real eignevalues are computed by using the inverse interation method of 

j Wielandt (as discussed in Ref. l8). 
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NUMERICAL EXAMPLES 

The results of two free vibration problems are presented to give an 
indication of the accuracy of the present method of analysis and to show the 
effects of the variation of certain parameters. The numerical examples con- 
sidered are a typical propeller blade euid a tapered beam having lateral dis- 
placements in only one direction. To verify the applicability of the equations 
of motion to a practical problem, the natural vibration frequencies of the 
propeller blade were determined numerically imd compared with experimental 
data. The computed results for the tapered beam are compared with exact and 
approximate solutions. 

Tapered Beams 

The free vibration characteristics of linearly tapered beams having 
lateral displacements in one direction are presented to substantiate the 
acciuracy of the solution method. Table 1 gives the comparison of computed 
and exact-solux ion bending modes of a nonrotating uniform beam at zero 
pitch. Table 2 gives a comparison cf bending frequencies of linearly 
tapered beams at zero pitch. Ihe nonrotating uniform beam results correspond 
to the exact values from Reference 19* The approximate Rayleigh-Southwell 
method presented in Reference 20 is used to estimate frequencies of the 
linearly tapered beaias. The computed frequencies of the linearly tapered 
beams are in excellent agreement with the Rayleigh-Southwell method. 
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j Propeller Blade 

^ The propeller blade selected for analysis is the WADC S-5 scale model 

of Reference 21. This blade was chosen since this reference gives a struc- 
tural description sufficient for the numerical solution as well as experi- 
. mental data for the natural vibration frequencies. The bl '.de is in effect 

cantilevered at O.IOI 6 meters from the center of rotation and the tip of 
the blade is at a radius of O .6096 meters. In the experimental program, 
tests were conducted for various pitch angles. The pitch settings were 
° defined by the values of 0 as measured at x = 0.75R where R is the 

I radius from the center of rotation to the tip of the blade. 

In order to compare numerical results with test data, solutions were 
computed for cases corresponding to the pitch settings and rotational speeds 

t • 

I of the experimental investigation. Numerical solutions were obtained by 

using 11 stations, which correspond to ten O.O 508 meter intervals, to describe 
the cross-sectional properties of the propeller blade. Physical properties 
, of the WADC S-5 blade, as given in Reference 22., are presented in Table 3. 

•Additional sectional properties needed for this analysis were estimated 
by assuming an elliptical cross section. These estimated physical properties 
are presented in Table k. 

The experimentally and analytically determined free vibration frequen- 
cies are given j.n Figures 2 and 3. Figure 2 illustrated the comparison for 
the first and second bending frequencies. Figure 3 shows the coniparison for 
the first torsion frequency. The computed results of Figure 3 were obtained 
by using the nonlinear twist euialysis of Appendix B. Fig:u'p 1* illustrates 

p 

B 

the effect of steady-state twist on the first torsion natural frequency. 

i 


Experimental torsion frequencies are compared with the nonlinear analysis of 
Appendix B and the linear analysis giv'en by the uncoupled fom of Equation 
(ic). Table 5 gives the percent error associated with each analysis. A 
comparison of the modal displacements is not possible since the mode shapes 
were not determined in the experimental investigation. 
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CONCLUDING REJ1APKS 

A numerical method for determi.iing the free vibration characteristics 
of a rotating cantilevered blade having nonuniform spanwise properties is 
presented. 3y expressing the equations of moti-jn in matrix notation, 
utilizing the integrating matrix as an operator, and applying the boundary 
conditions, tne equations are formulated into an eigenvalue problem vhose 
solutions may then be determined by conventional methods Computed results 
were compared with experimental and analytical data. The comparison indicates 
that the method of solution yields very accurate results. 
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APPENDIX A 
EQUATIONS OF MOTION 

The eque.tions of motion for combined bending and torsion are obtained 
by application of the theory developed by Houbolt and Brooks.^ This case 
represents the coupled bending and torsion of a tirtsted rotating beam where 
the elastic axis, mass axis, and tension axis are not necessarily coincident. 
Figure 5 depicts the position of the elastic axis as a space ciurve. The 
applied loadings are shown acting at a radied distance (n) from the axis 
of rotation. The bending moments produced by these loadings at a radial 
position (x) are 


M 


/{- 


p^[w(n) - 


w(x)] + 


pjn - 


:] + 




dn 


(Ala) 


M 

z 




y - x] + dn (Alb) 

/• R 

J { lx " lx ^ " Py(('^(n) - '^(*)) 


- (n - x) 1^ ] + P^[(v(ti) - v(x)) - (n - x) ] 


(Ale) 


Mqx + qy 


9v 

3x 


'L 


3w J 

3x f 


dn 
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For the present analysis, blade elasticity may be adequately described 
by the conventional bending end torsion characteristics developed in 
Reference 1. Thus, the elastic restoring moments at an arbitrary radial 
position (x) are 


M * (EI^co 8^0 + EI-sin^6)rf'' + (EI^ - El, )sin 0 cos 0 v" 
y 1 2 2 1 


(A2a) 


I 1 

- (Te + EB 0 <p ) sin 0 - Te 4> cos 0 
a 2 a 


M = (El, sin^0 + EI„ cos^ 0)v" + (EI^ - El ) sin 0 cos 0 W 

« f I 

- (Te + EB- 0 ) cos 0 + Te sin 0 (A2b) 

a 2 a 

M. = [GJ + Tk ^ + EB.(0')^] 4)' + Tk ^ 0* 

U & X & 

(A2c) 


- EBg 0 (v" cos 0 + w" sin 0) 


Equation (Ale ) may be written as 


3M 

= ,.q - q. V - q w* - M V*' + M w" 

3x X y z y z 


(A3) 


substituting Equations (A2a) and (A2b) eind neglecting products of elastic 
deflections yields 


3M 

t _ 

3x Sc 



q w' + Te v" sin 0 - Te w" cos 0 
z a a 


(AM 
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Combining Equations (Al), (A2), (aU) and substituting the applied loadings 
given by Equation (21) of Reference 1 yields 


(EI^ cos^ 6 + Elg sin^ 0)w" + (EI^ - EI^) cos 0 sin 0 v" 


- EB- 0'((»* sin 0 - Te ^ cos 0 

b £t 




{mn[w(n) - w(x)] 


/•R 

+ men it* cos 0} dn = - / (mi? + me$ cos 0) (n - x)dn 


(Elg - EIj^)cos 0 sin 0 w" + (EI^ sin“^ 0 + cos 0)v" - EB^ 0'4>* cos 0 


+ Te 4 sin 0 + n 
a 


r <4 


xv(n) - nv(x)] - men<J) sin 0} dq 


(A5b) 


-r 


(mi? - me$ sin 0) (q - x)dq 


[GJ + Tl; ^ + EB- (0')^](<)' - EB- 0' (v" cos 0 + w" sin Q) + f Te (v" sin 0 
a 1 d Jo. 

X 

2 ^ 2 2 

- w" cos 0) dq + / [mev sin 0 + m (k - k )<*> cos 20 

•'x “2 “l 

- meq (v' sin 0 - w' cos 0)] dq = - / [me(i? cos 0 - i? sin 0) 


+ 4*m k ] dq 
in 
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APPENDIX B 

STEADY-STATE TWIST OF ROTATING BLADES 

Equation (l, is valid for small displacements from the undeformed configu- 
ration of the nonrotating blade. The angle 0 at any radial station x, 
defines the orientation of the Jiajor principal axis relative to the plane of 
rotation. Thus, the local orientation of the nonrotating blade includes 
pitch and pretwist. However, the local orientation of a blade segment 
may be modified due to the presence of rotation. This effect is defined as 
"steady-state" twist and is a function of blade characteristics and rotational 
speed. 

An approximate estimate of the steady-state twist may be obtained by 
assuming that blade torsional response is uncoupled. The uncoupled linear 
torsional restoring moment is given by 

M, = [GJ + Tk^ + EB, (6' + Tk^ 9' (Hi) 

Equation (Bl) is valid for v = w = 0, and 4> a small angle. However, the 
presence of steady-state twist may invalidate the assumption of small angular 
displacement from the undeformed position. The equivalent form of Equation 
(Bl) for V = w = 0, and large angular displacement is (see Appendix A, 

Ref. 1) 

M. = {GJ + EBj(6')^ + I 6'<J" + ^ (<t’ )^]}<J>' + Tk^ O' + <J>') (B2) 

where $ is the sum of the steady-state twist, (1)^, and the torsional elastic 




deformation ({>. 
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4>(x,t) « * 4>(x,t) (B3) 

Using Equation (B2) the uncoupled torsional equation is 

{GJ + EB-i(6' ^ 0'4>' + J ♦ Tk^ (6' ♦4'') 

X 2 c & 

mk^ 4> dn = 0 (BU) 

Cl 

Substituting ‘Equation (B3) into Equation (BU) and setting all pei^'.urbation 
quantities equal to zero yields the equilibrium equation 




m(k_^ - k_^) sin 2(0 + 4>) dn + J" 




m. 


{GJ + EBJ(0')^ + I 0'<l>' + i (4>')^]}4>: + Tkf (0* + (^') 

1 2o2ooa o 

^ B p p 

+ ^ / m(k - k ) sin 2(0 + 4 ) dn = 0 (B5) 

2 y. m- m o 


Substituting Equation (B3) into Equation (BU), subtracting the equilibrium 
equation and discarding higher order products of perturbation quantities, 
yields the perturbation equation. 

{GJ + Tk^ + EBJ(0')^ + 30'4>' + I (<|)')^]}<^' 
a 1 o 2 o 

22 /*^ 2 ’* 

+ Q I m(k - k )4> cos 2(0 + (|))dn+ / mk ({idn = 0 (B6) 

J m_ m, o y m 

X 2 1 X 


The steady-state twist (|) is determined by iteration from Equation (B5). 

o 

This value is used to determine the coefficients of Equation (B6) which Is 
solved for the uncoupled torsional frequencies. 
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TABLE 2. BENDING FREQUENCIES OF LINEARLY TAPERED BLADES 
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TABLE 3. PHYSICAL PROPERTIES OF PROPELLER BLADE 







TABLE U. ESTIMATED PHYSICAL PROPERTIES OF PROPELLER BLADE ASSUMING ELLIPT'ICAL 













TABLE 5. CO^^PARISON OF EXPERD-IENTAL AND ANALYTICAL TORSION FREQUENCY 
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Figure 2. - Variation of bending frequencies with rotational speed. 




Figure 3. Variation of torsion frequency with rotational speed 



Figure 4. - Torsional frequ-'ncy versus rotational speed 
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Figure U. - Concluded 




Figure 5.- Equilibrium of forces and moments. 





